function [ gN ] = gravity_qyy( CNE, h, RE ) %#codegen
%GRAVITY_QYY 按照秦永元老师课本模型计算重力加速度
%
% Tips:
% # E系EUGw，初始N系ENU（可以为非地理坐标系）
%
% Input Arguments:
% # CNE: 3*3矩阵，位置矩阵
% # h: 标量，高度，单位m
% # RE: 标量，地球长半轴，单位m
%
% Output Arguments:
% # gN: 3*1向量，N系下的重力加速度，单位m/s^2
%
% References:
% # 秦永元. 惯性导航[M]. 科学出版社. 2006
%

sinL = CNE(2, 3); 
sinL2 = sinL*sinL;
gt = 9.7803267714 * (1 + 0.00193185138639*sinL2) / sqrt(1 - 0.00669437999013*sinL2) * RE*RE/(RE+h)/(RE+h); % REF1式7.2.17 - 7.2.18
gN = [0; 0; -gt];
end